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ABSTRACT 

We present the first ever direct TV-body computations of an old Milky Way globular 
cluster over its entire life time on a star-by-star basis. Using recent GPU hardware at 
Bonn University, we have performed a comprehensive set of TV-body calculations to 
model the distant outer halo globular cluster Palomar 14 (Pal 14). Pal 14 is unusual in 
that its mean density is about ten times smaller than that in the solar neighborhood. 
Its large radius as well as its low-mass make it possible to simulate Pal 14 on a star-by- 
star basis. By varying the initial conditions we aim at finding an initial TV-body model 
which reproduces the observational data best in terms of its basic parameters, i.e. 
half-light radius, mass and velocity dispersion. We furthermore focus on reproducing 
the stellar mass function slope of Pal 14 which was found to be significantly shallower 
than in most globular clusters. 

While some of our models can reproduce Pal 14's basic parameters reasonably 
well, we find that dynamical mass segregation alone cannot explain the mass function 
slope of Pal 14 when starting from the canonical Kroupa initial mass function (IMF). 
In order to seek for an explanation for this discrepancy, we compute additional initial 
models with varying degrees of primordial mass segregation as well as with a flattened 
IMF. The necessary degree of primordial mass segregation turns out to be very high, 
though, such that we prefer the latter hypothesis which we discuss in detail. This 
modelling has shown that the initial conditions of Pal 14 after gas expulsion must have 
been a half-mass radius of about 20 pc, a mass of about 50000 M , and possibly some 
mass segregation or an already established non-canonical IMF depleted in low-mass 
stars. Such conditions might be obtained by a violent early gas-expulsion phase from 
an embedded cluster born with mass segregation. Only at large Galactocentric radii 
are clusters likely to survive as bound entities the destructive gas-expulsion process 
we seem to have uncovered for Pal 14. 

In addition we compute a model with a 5% primordial binary fraction to test if 
such a population has an effect on the cluster's evolution. We see no significant effect, 
though, and moreover find that the binary fraction of Pal 14 stays almost the same 
and gives the final fraction over its entire life time due to the cluster's extremely low 
density. Low-density, halo globular clusters might therefore be good targets to test 
primordial binary fractions of globular clusters. 

Key words: globular clusters: Palomar 14 - initial mass function - methods: N-body 
simulations 
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1 INTRODUCTION 

The initial mass function (IMF) of stars is an important 
quantity for astrophysics. Since the properties and evolu- 
tion of a star are closely related to its mass, the IMF is 
a link between stellar evolution and the evolution of stel- 
lar systems. It influences the chemical evolution of galaxies 
and the dynamical evolution of star clusters. The IMF is 
often expressed as a power-law function (4^- oc m~ a ). The 
traditional IMF that was proposed by Salpeter (1955) from 
analyzing stars in the solar neighborhood with masses be- 
tween O.4M0 and IOMq, has a = 2.3, implying that high 
mass stars are rare compared to low mass stars. Later ad- 
vances showed that the IMF slope flattens to a = 1.3 for 
stars in the mass range 0.5 Mp) — 0.08Mq while it remains 
Salpeter's for massive stars (|Kroupall200lh . 

Since stars in star clusters have similar metallicities, 
distances and ages, they are the best objects to study the 
stellar mass function. However, the mass function of stars in 
clusters evolves through stellar and dynamical evolution and 
one has to understand how these processes have affected the 
mass function of stars before one can extract the IMF from 
the observed mass function (Baumgardt & Makino 2003). 
For example, as a result of energy equipartition, the heavier 
stars tend to move toward the cluster centre while lighter 
stars tend to move further away from the cluster centre. 
Since observations of stars in globular clusters usually are 
radially limited, the effect of mass segregation has to be con- 
sidered when deducing the global mass function from the 
observed local one. Moreover, as a result of preferential loss 
of low-mass stars, a cluster's mass-to-light ratio evolves dif- 
ferently from that expected from pure stellar evolution, such 
that dyna mical evolution leads to a decrease of the mass-to- 
light ratio ( Baumga rdt fc Makinoll2003l ; iKruiissen fc Mieskei 
|2009| ; lAnders et alj|2009h . 

While the observed mass function shows local varia- 
tions, e.g. as a consequence of dynamical mass segregation, 
the IMF is usually assumed to have no spatial variation 
throughout a cluster. However star formation simulations 
indicate that this assumption could be invalid (e.g., Tan et 
al 2006; Krumholz et al 2009). There is also some observa- 
tional evidence for deviations from the standard initial mass 
function, at low and high mass ends, in many star clusters, 
which might be the result of dynamical evolution (e.g., see 
Elmegreen 2004 for review). 

Some of the ideas that have been proposed to explain a 
shallowness of the slope at the high mass end include primor- 
dial mass segregation of stars in the cluster (e.g., Vesperini 
& Heggie 1997; Kroupa 2002; Mouri & Taniguchi 2002). But 
so far it is not clear if star clusters were primordially mass 
segregated from the star formation epoch or obtained mass 
segregation via dynamical evolution at a later time. Espe- 
cially for some young star clusters it has been shown that the 
observed mass segregation cannot be explained dynamically 
from an initially unsegregated system (Bonnell & Davies 
1998). The early mass loss due to stellar evolution in star 
clusters leads to shorter lifetimes, and faster expansion, and 
has a stronger impact on initially segregated clusters than 
on unsegregated clusters. (Vesperini et al. 2009b). 
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Measured mass functions of a sample of globular clus- 
ters surprisingly show that all high concentration clusters 
have steep mass functions (i.e., larger a), while low concen- 
tration ones have a smaller a (De Marchi et al. 2007). This 
is puzzling, because concentrated clusters are expected to be 
dynamically more evolved and thus most of their low mass 
stars should have evaporated from the clusters. However, 
Marks et al. (2008) showed that for initially mass-segregated 
clusters, mostly low-mass stars are lost due to gas expulsion 
which implies a shallower slope in the low-mass range. In- 
deed this approach allows to deduce the evolution of the 
very early Milky Way in unprecedented time resolution as it 
allows us to deduce the changes in the cluster potential in- 
duced by gas-expulsion and the varying Galactic tidal field 
l|Marks fc Kroupall2010t ). 

Several studies have been done to dynamically model 
the evolution of star clusters using iV-body simulations. Re- 
cently, Harfst et al. (2009) modelled the very young Arches 
cluster to find the best fitting initial model by comparing 
simulations with observational data. In this way they aimed 
at constraining the parameters for the initial conditions of 
the cluster such as IMF, size, and mass of the cluster by ne- 
glecting the Galactic potential and stellar evolution. More- 
over, using iV-body simulations, a close correlation between 
the slope of the low-mass end of the stellar mass function and 
the fraction of the initial cluster mass loss has been found 
by several authors (Vesperini & Heggie 1997; Baumgardt & 
Makino 2003; Trenti et al. 2010). 

In the present paper we report on simulations of the 
diffuse outer halo globular cluster Palomar 14 (Pal 14). Pal 
14 with a mass of about 10 4 M© and a (three dimensional) 
half-light radius of 34 pc, is a particularly interesting object 
to study. It is one of the most extended Galactic globular 
clusters and, with a distance of ~ 70 kpc (Jordi et al. 2009) 
is also among the outermost clusters of the Milky Way. This 
large Galactocentric distance and its low stellar mass density 
make this cluster an ideal candidate to decide between New- 
tonian gravity and MOND (Baumgardt et al. 2005, Haghi 
et al. 2009). Its low density suggests that binary stars may 
be an impo rtant issue for interpretin g its measured velocity 
dispersion (|Kiipper fc Kroupall201(il ). Its low mass together 
with its large radius make it also possible to simulate Pal 14 
on a star-by-star basis. However, simulations of the evolu- 
tion of individual globular cluster's long-term evolution have 
so far been done only with Monte Carlo or Fokker Planck 
simulations l|Giersz fc He ggie 2009) . We therefore performed 
a set of direct iV-body simulations to determine the most 
likely starting conditions of Pal 14 in terms of total mass, 
initial half-mass radius and stellar mass function. We also 
test if primordial mass segregation is necessary to reproduce 
the observed flattened mass function. Moreover, we check if 
primordial binaries get frequently disrupted within such a 
sparse cluster or whether Pal 14 has preserved its primor- 
dial binary content. It should be noted that all simulations 
performed in the present study start once the cluster is gas- 
free and has settled back into virial equilibrium. 

This paper is organized as follows: in Section 2 we de- 
scribe the observational data we use to compare our simula- 
tions to. Our methodology is described in Section 3. In Sec- 
tion 4 we explain the different N-body set-ups and compare 
the numerical results with the observed data. Conclusions 
are given in Section 5. 
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2 OBSERVATIONAL DATA 

We use observational data of Hilker (2006) and Jordi et al. 
(2009) who have presented a spectroscopic and photometric 
study of Pal 14. Jordi et al. (2009) have obtained the clus- 
ter's mass function using HST/Wide Field Planetary Cam- 
era 2 (HST/WFPC2) data from the HST archive. Moreover, 
they determined Pal 14's velocity dispersion by measuring 
the radial velocities of 17 giant stars using the Ultraviolet- 
Visual Echelle Spectograph (UVES; Dekker et al. 2000) at 
the Very Large Telescope (VLT) of the European Southern 
Observatory (ESO) in Chile and the High Resolution Echelle 
Spectrograph (HIRES; Vogt et al. 1994) on the Keck I tele- 
scope. 

The distance of Pal 14 from the Sun was derived from 
the apparent magnitude of the horizontal branch in the 
colour-magnitude diagram and an isochrone fit to be about 
71 ± 1.3 kpc (Jordi et al. 2009). This is closer to the Sun 
than previous estimates, e.g. Hilker (2006) and Dotter et al. 
(2008) derived larger distances of 74.7 kpc and 79 kpc, re- 
spectively. Pal 14 has a projected half-light radius of l'.28" 
(Hilker 2006), which was derived from a curve-of-growth 
analysis of the integrated light of the member stars. This cor- 
responds to a projected half-light radius of R p hi = 26.4 ±0.5 
pc at the assumed distance of 71 ± 1.3 kpc and implies a 3D 
half-light radius of about 35.6 ± 0.6 pc. 

Jordi et al. derived the age of the cluster to be 11. 5± 0.5 
Gyr by adopting [Fe/H] = -1.50 for the metallicity (Harris 
1996) from the best isochrone fitting of Dotter et al. (2007). 
Although previous work has shown that Pal 14 is 3-4 Gyr 
younger than typical halo GCs with the same metallicity, 
thus has an age of ~ 10 Gyr (Sarajedini 1997; Hilker 2006; 
Dotter et al. 2008), the recent estimated age by Jordi et al. 
reduces this difference. 

The density profile of the cluster does not exhibit a 
significant ellipticity and the distribution of its stars is com- 
patible with a homogeneous, circular distribution within the 
uncertainties (Hilker 2006). The core radius, and the central 
surface brightness of the cluster were determined from the 
best fitting King (King 1966) profiles. The concentration 
parameter which is defined as c = log(r t /r c ) is about 0.85 
(Hilker 2006), i.e. rather low compared to most Milky Way 
globular clusters. 

The velocity dispersion of Pal 14 by measuring radial 
velocities of 17 red giant cluster members has been presented 
with and without including one probably binary or AGB 
star (hereafter star 15), whose velocity is more than 3<r 
away from the mean of the other stars (Jordi et al. 2009). 
For the sample including star 15, the global line-of-sight 
velocity dispersion is 0.64 ± 0.15 kms -1 , and without star 
15 it is 0.38 ± 0.12 kms" 1 (Jordi et al. 2009). The quoted 
uncertainties reflect the 1-er error. 

The upper boundary of Pal 14's main sequence corre- 
sponds to a magnitude of 23.4 mag at a wavelength of 555 
nm. Based on photometric data, Jordi et al. (2009) used 
the masses given by the 11.5 Gyr isochrone by Dotter et al. 
(2007), and found that the observed stellar masses on the 
MS covered the mass range of 0.495 M© up to 0.825 M . 
The mass function of Pal 14 was found to be flatter than 
a Kroupa or Salpeter IMF in the observed mass range, and 
the best fitting power-law mass function was obtained for a 
slope of a = 1.27 ± 0.44. 




Figure 1. Colour-magnitude diagram of Pal 14. The grey dots 
are all photometric data, the black dots are a selection of point 
sources with reasonable photometric errors (< 0.2 mag). The lines 
give the selection regions of RGB and HB stars. We put the bot- 
tom limit of the RGB to Vb = 22.9 mag, which would correspond 
to Mvfi = 3.48. Red squares are those stars in the RGB region 
that are within the tidal radius of Pal 14. The same applies for 
the blue squares for the HB stars. 

The total mass of stars inside the half-light radius and 
within the mass range 0.525 M to 0.825 M has been esti- 
mated to be M = (2200 ± 90)M Q . With the measured slope 
for masses down to 0.5 M and using a Kroupa-like mass 
function, a — 1.3 for masses between 0.1 M and 0.5 M , 
Jordi et al. (2009) estimated the total mass of Pal 14 to be 
about Mtot = 12OOOM without taking into account stel- 
lar remnants. Since this estimated total mass is a derived 
quantity and not directly measured, we use the number of 
bright main-sequence stars, Nt, s to compare with our mod- 
els, instead. For this aim, we used the star-count data with 
completeness factor greater than 0.50, which leads to the 
mass range from 0.525 M© to 0.795 M© and consequently 
to the number of Pal 14's bright stars N bs = 2954 (Table 
3 of Jordi et al. (2009)). This mass range corresponds to a 
magnitude range 27.2 > 771555 > 23.5. 771555 is the HST filter, 
F555W, that is very similar to Johnson V (Fig. 7 of Jordi et 
al. (2009)). The stellar mass at the main-sequ ence turn-off 
point , rriMSTOfi = 23.44, is about 0.79 M (|Dotter et alj 
120071 ). Consequently the central and average density of the 
cluster are 0.1 and 0.01 stars/pc 3 , respectively. The two- 
body relaxation time of Pal 14 is about 50 Gyr, i.e. several 
times larger than a Hubble time. 

The total number of giant stars including the red gi- 
ant branch (RGB) and horizontal branch (HB) stars is also 
another observable quantity. The photometry published by 
Hilker (2006) has not covered entirely the cluster (a segment 
in the West was missing) and was not corrected for fore- 
ground/background contamination. Therefore, the numbers 
in the tables of Hilker (2006) are not complete. We there- 
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Figure 2. Spatial distribution of the photometrically observed 
Pal 14 stars. The big rectangle marks the area of the observations. 
The red dots are all objects in the RGB regime shown in the CMD, 
the blue dots are the objects in the HB area of the CMD. The 
two large circles are the half-light radius at 1.15 arcmin and the 
tidal radius of 4.81 arcmin (both from Baumgardt ct al. (2005)). 
The region is divided into stripes/segments, called EE, E, W and 
WW from East to West. Everything outside the tidal radius is 
the background (B) which was used to count the contaminating 
objects. 

fore redo a new selection and incompleteness calculation to 
count all RGB and HB stars. Figure [1] shows Pal 14's colour- 
magnitude diagram (CMD). Note that some of the brighter 
RGB stars on the bluer side might actually be AGB stars. 
This cannot be differentiated from this photometry. From 
this CMD it also becomes clear that it makes no sense to 
count SGB stars. The photometric error is large at those 
magnitudes and the counts are highly incomplete. The spa- 
tial distribution of the observed stars is shown in Figure [2] 
We divided the regions into stripes/segments, called EE, E, 
W and WW from East to West. Everything outside the tidal 
radius is the background (B) which was used to count the 
contaminating objects. Obviously, the missing region WW 
should contain more or less the same number of RGB/HB 
stars as the region EE. But one can also make other sanity 
checks of number counts, for example whether there are the 
same amount of RGB/HB stars in the E and W stripes, or 
whether really the half-light radius divides the number of 
RGB/HB stars into half (assuming no mass segregation). 

The total area within the tidal radius (r < r t ) is 78.5 
arcmin 2 . There are different ways to derive the number of 
RGB and HB stars: 1) assuming that there are as many RGB 
stars in the WW segment as in the EE segment, 2) doubling 
the number of RGB stars in the East (E+EE). The first 
method gives Nrgb = 164T18, the second Nrgb = 169±24, 
agreeing very well with the first number. Both values give 
the background subtracted number counts, taking into ac- 
count the ratio of the selected area vs. the background area, 



i.e. scaling the background counts to the selected area. This 
analysis shows that one can safely assume that the RGB 
counts in the West should be similar to those in the East. 
Another result is that there are slightly more RGB stars in- 
side the half-light radius, Nrgb = 46 ± 7, than outside the 
half-light radius, Nrgb = 39 ± 10. The numbers are consis- 
tent within the errors, but might hint to mass segregation 
of RGB stars on the East side. By the same analysis for 
the HB stars in the tidal radius we find Nhb = 34 ± 7 for 
method 1) and Nhs = 33 ± 10 for method 2), both val- 
ues being in very well agreement. The total number of giant 
stars (RGB+HB), N 9 = 198 ± 19, is given in the tables to 
compare with those obtained from the modelled clusters. 

In order to reproduce those observational values, e.g. 
the half-light radius, R p hi , the number of bright stars, jVj, s , 
and the mass function slope, a, we construct a set of JV-body 
models for Pal 14 in the next section. 



3 JV-BODY MODELS 

In total we computed a set of 66 models to find the ini- 
tial conditions which reproduce the observations of Pal 14 
best. The simulations were carried out with the collisional 
jV-body code NBODY6 (Aarseth 2003) on the GPU com- 
puters of the Stellar Population and Dynamics Group of the 
University of Bonn. 

All clusters were set up using the publicly available code 
McLuSTEEQ(Kiipper et al. (2010), in preparation). We used 
a Plummer profile in virial equilibrium as initial mass dis- 
tribution and computed all clusters for 11 Gyr. We chose 
the initial mass and size (i.e. half-mass radius) of the clus- 
ter as free parameters and kept the other parameters of the 
initial conditions fixed. The number of stars was varied in 
the range 7 x 10 4 < N < 10 5 such that the clusters have 
total masses in the range of 4.6 — 5.4 x 10 4 Mq. Since the 
half-mass radius grows due to stellar-evolutionary mass loss 
and dynamical evolution over the cluster life time, we chose 
the initial 3D half-mass radii in the range of 18 pc to 22 pc 
to reach projected half-light radii of about 26 pc after 11 
Gyr. 

We have chosen two different initial mass functions for 
our models. First, we start with the canonical IMF (Kroupa 
2001) using lower and upper mass limits of 0.08 Mq and 
100 Mq, respectively. At the second step, the clusters were 
constructed using a flatter mass function with lower and 
upper mass limits of 0.08 Mq and 5 Mq, respectively, mim- 
icking the initial conditions of an initially mass-segregated 
star cluster after gas expulsion at an age of about 100 Myr. 
Stellar evolution was modelled according to the routines by 
Hurley et al. (2000). Due to the low escape velocity of Pal 
14, we assume a 0% retention fraction for neutron stars and 
black holes which form during the simulation. Note that the 
canonical (universal) two-part power-law IMF has ai = 1.3 
for 0.08 < m < 0.5 M and a 2 = 2.3 for m > 0.5 Mq. Here 
we refer to a as the index of the mass function in the stellar 
mass range 0.525 M© to 0.795 MQ. 

For simplicity, and since no information on the proper 
motion is available, the clusters move on a circular orbit 
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through a logarithmic potential, <f>(Ra) — Vq InRc, of the 
Milky Way at the assumed distance of Pal 14 and with a 
circular velocity Va = 220 km/sec. 

All but one cluster do not have primordial binaries. We 
constructed one model with a primordial binary fraction of 
5% to study the impact of binaries on cluster evolution and 
to estimate the fraction of binaries which dissolve in a Hub- 
ble time in a cluster as sparse as Pal 14. The binary popu- 
lation was set up using the Kroupa (1995) period distribu- 
tion with random pairing of the binary companions. Ran- 
dom pairing of companion masses for primaries less massive 
than 1 Mq is the correct procedure to reproduce the pre- 
main sequence and main sequence binary-star data (Kroupa 
2008). 

The modelled clusters were assumed to be not initially 
mass segregated. In addition we computed a number of clus- 
ters with different degrees of primordial mass segregation to 
investigate the effect on the cluster's evolution. An overview 
of the clusters is given in Tables 1, 2 and 3. 



4 FINDING THE BEST FITTING MODEL 

In order to find the most likely initial conditions, a grid 
of models with different initial half-mass radii and initial 
masses was calculated. The main observational values which 
were used to compare with each model, are the present day 
number of bright main-sequence stars, Nbs, the present day 
projected half-light radius, R p hi, and the slope of the mass 
function in the mass range 0.525 M© to 0.795 M© within 
the half-light radius, a, with values 

N bs = 2954 ±175, (1) 
Rphi = (26.4 ± 0.5) pc (2) 
a = 1.27 ±0.44, (3) 

Note that, to match R p hi in the iV-body models we 
use only the giant stars since the projected half-light ra- 
dius comes from Hilker (2006) and he could only see the 
giant stars in Pal 14. Moreover, for iV;, s we count the main- 
sequence stars inside the half-light radius with masses be- 
tween 0.525 M Q and 0.795 M Q . Finally, we extract the 
mass function for those main sequence stars within the half- 
light radius, with a mass in the range of 0.525Mq < m < 
O.795M0, i.e. in the same way as Jordi et al. have done. We 
apply 10 bins between log 10 m = —0.1 and log 10 m = —0.3. 
We reject the lowermost point in the observational data of 
Jordi et al. since it has a relatively high incompleteness cor- 
rection, i.e. uncertainty. Since compact remnants are unob- 
servable in Pal 14, as they are too faint, we ignored them in 
our measurement of the mass function. 

To compare our models in an objective way and to mea- 
sure the quality of the fit of the model to the observations, 
we define a goodness-of-fit parameter, g p . This parameter 
depends on the model parameters, and measures the close- 
ness between the observations and the model as follow 

i Pobservation Pmodel / .\ 

9p = 1 , (4) 

P 'observation 

where p represents the corresponding model parameter. The 
goodness-of-fit parameter ranges from to 1, where 1 de- 
scribes an excellent fit. Since we have three main observ- 



ables, Rphi, Nbs and a, we also have three goodness-of-fit 
parameters, ga, gjv and g a , respectively. 

There is a significant statistical dispersion in the re- 
sults if we repeat the calculations for one arbitrary model 
with the same initial parameters but with a different ran- 
dom seed. The required computational GPU time doesn't 
allow us to do statistics for all parameters, though. Hence, 
in order to estimate the errors of the resulting parameters, 
we do several runs (10 times) for a particular model and dis- 
cuss the resulting uncertainty. In the following investigation 
we will use these uncertainties for our discussion. 

In the following sections we will investigate four sets of 
models. First, we investigate the initial conditions which we 
call 'regular'. We assume a Kroupa IMF without primordial 
mass segregation or primordial binaries. In the second step 
we assume that the IMF was flattened within the first 100 
Myr through stellar evolution and gas expulsion processes. 
The third set of models is computed to test how much pri- 
mordial mass segregation we have to add in order to explain 
the observed mass function of Pal 14. In the last section 
we test the influence of primordial binaries on the cluster 
evolution and on the evolution of the mass function. 

4.1 Regular initial conditions 

Figure [3] shows how the goodness-of-fit parameters of Nbs, 
Rphi, and a vary for the different values of the model pa- 
rameters, initial half-light radius and initial cluster mass. 
Moreover, the product of the goodness-of-fit parameters for 
N os and R p hi is shown in the lower right panel. 

The darker areas show where the goodness-of-fit pa- 
rameters are close to 1 in the parameter space, i.e. show the 
range in initial conditions that best fit to Pal 14. We find 
that many models (e.g. M54R22, M52R21, and M50R20) 
can reproduce Pal 14's half-light radius and its number of 
bright main-sequence stars. 

We look for models whose velocity dispersion lies within 
the uncertainties of the observational values. Therefore, we 
take only stars within a clusters tidal radius and with a mass 
higher than 0.8 Mq, and measure their velocity dispersio n 
and compare it to the observed value of Ijordi et al.l (2009). 
According to Table [T] the line-of-sight velocity dispersion in 
most of our models is about 0.6 km s - . 

Given the small sample size (i.e., radial velocity of 17 
stars are used to determine the observed velocity disper- 
sion), the observed velocity dispersion has a large error 
which could be much larger than the formal error (Gentile 
et al. 2010). In order to see how well a velocity dispersion 
of ~ 0.6 km s~ agrees with the observed one, we use the 
Kolmogorov-Smirnov (KS) test, which is the applicable sta- 
tistical test for such a small sample size. 

In Fig. [4] we show the cumulative distribution function 
(cdf) vs. radial velocity for both samples with and with- 
out star 15, separately. This cdf we compare to the cdf of 
one of the best fitting models (M50R20). With star 15, the 
maximum difference between the two cdfs is 0.13, which cor- 
responds to P-value of 0.94, while in case without star 15 the 
maximum difference between the two cdfs and correspond- 
ing P-value are 0.14 and 0.91, respectively. This means that 
the models can be rejected with 6% (including star 15) and 
9% (excluding star 15) confidence. 

In addition, we find that the mean value of the velocity 
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Figure 3. Goodness-of-fit parameters from the measured half-light radius (upper left panel), number of bright stars (upper right panel) 
and mass function slope (lower left panel) for the models with regular initial conditions. In the lower right panel the product of the 
two upper panels is given. The varying model parameters are the initial half-mass radius and the cluster initial mass. A goodness-of-fit 
parameter close or equal to 1 is best and indicated by dark shaded areas in the plots. 



dispersion of all stars is similar to that of the red giant stars 
and main sequence stars with masses larger than 0.8 Mq. 
Therefore, the velocities of the giants are representative of 
all stars in the cluster. The number of giant stars, N g , is also 
in most of the modelled clusters in good agreement with the 
observations. 

Figure [5] shows the mass function of one example model 
before and after the computation. The red and green points 
denote the number of stars inside and outside the half-light 
radius, respectively. We see that after 11 Gyr the mass func- 
tion inside the half-light radius is significantly different from 
the one outside the half-light radius, which implies that mass 
segregation has taken place. The mass function of the stars 
inside the half-light radius is clearly flatter than those of the 
stars outside. It seems that even in a cluster as large as Pal 
14, dynamical mass segregation plays an important role. 

According to Table [1] and lower left panel of Figure [3] 
the slope of the mass function, a, in all models is far from the 
observed value of Pal 14, though, even when accounting for 
observational errors. The upper limit of the observed value 
within the 1-a error is a — 1.71, and only if we account for 
the statistical error for the modelled slope, which is ~0.15, 
then two models (i.e., M48R19 and M50R18) are marginally 
compatible with the observed slope. 

Our simulations have typically a ~ 2.0 after 11 Gyr , 



i.e. only a mild decrease of the low-mass-star slope from the 
initial Salpeter value of a = 2.35 due to mass segregation. In 
contrast, the observed slope is a = 1.27±0.44, which agrees 
with this value only at the 1.75 a level. There is therefore 
some evidence that Pal 14 started with a mass segregated 
Kroupa IMF. 

In order to see whether an extension of the initial pa- 
rameter space would resolve the problem, we repeat the sim- 
ulation for some extremely high and low initial values, i.e. we 
chose 25 pc and 15 pc for the half-light radius and 40000 Mq 
and 60000 Mq for the initial mass. The results are shown 
in Table Q] As can be seen, even by selecting such values 
for the initial conditions, the present day mass function of 
Pal 14 cannot be reproduced. Moreover, we find that the 
values for jVj, s and R p hi get much worse in these directions 
of the parameter space, showing that we were searching in 
the right place initially. 



4.2 Flattened mass function 

As shown in Section [4.11 the initially unsegregated clusters 
with a canonical Kroupa IMF cannot reproduce the observed 
slope of the mass function well. The observed slope is a — 
1.27 ± 0.44 but the slope that could be achieved from N- 
body calculations is about 2. In this section we show that 
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Table 1. Details of all iV-body models with regular initial conditions (Sec. 14, ljl . The first column gives the 
model name, in which the first two digit numbers denote the initial mass in 1000 Mq and the second two digit 
numbers denote the initial half-mass radius in parsec. Columns 2 and 3 are the projected half-mass and half-light 
radius, Rphm an d R p hl> respectively. Column 4 is the present day number of bright main-sequence stars inside 
the half-light radius, with masses between 0.525 and 0.795Mq, N b3 , whereas column 5 gives the corresponding 
total mass of these stars, a is the present day slope of the mass function for stars with masses between 0.525 
and O.795M0. Column 7 is the line-of-sight velocity dispersion. The last column gives the total number of 
giant stars in the cluster including RGB and HB stars. Typical errors of the numerical values were obtained by 
repeating runs for a particular model and are indicated in the header. In the last line, the observational values 
are given for comparison. The subscript '*' indicates the observed velocity dispersion with taking star 15 into 
account. Note that none of these models are an acceptable fit to Pal 14. 



Model 


Rphm [pc] 


Rphl [pc] 




M r<R h J M ^ 
(±320) 


a: 


oi os [km/sec] 






(±0.8) 


(±2.3) 


(±260) 


(±0.15) 


(±0.01) 


(±16) 


M40R15 


19.6 


15.6 


2374 


1521 


2.10 


0.65 


162 


M40R25 


31.5 


30.3 


2182 


1402 


1.85 


0.51 


136 


M46R18 


23.5 


22.0 


2912 


1869 


1.96 


0.64 


163 


M46R19 


24.3 


19.8 


2382 


1530 


1.90 


0.64 


161 


M46R20 


25.8 


23.8 


2661 


1705 


2.04 


0.61 


174 


M46R21 


28.4 


25.2 


2407 


1542 


2.08 


0.58 


188 


M46R22 


29.3 


26.8 


2544 


1627 


2.10 


0.58 


160 


M48R18 


23.7 


23.1 


3311 


2127 


1.88 


0.66 


174 


M48R19 


24.3 


21.5 


2847 


1821 


1.77 


0.64 


160 


M48R20 


26.2 


22.5 


2641 


1689 


2.13 


0.61 


172 


M48R21 


28.0 


25.6 


2687 


1722 


2.02 


0.59 


174 


M48R22 


28.1 


27.6 


3079 


1970 


2.14 


0.61 


171 


M50R18 


24.2 


21.3 


2992 


1925 


1.78 


0.68 


194 


M50R19 


24.5 


19.5 


2638 


1695 


1.87 


0.65 


193 


M50R20 


26.6 


25.8 


3108 


1995 


1.97 


0.64 


170 


M50R21 


26.9 


24.1 


2796 


1788 


2.14 


0.61 


178 


M50R22 


29.3 


29.0 


2997 


1918 


2.10 


0.60 


163 


M52R18 


24.3 


21.5 


3181 


2038 


2.06 


0.68 


181 


M52R19 


25.7 


24.0 


3245 


2079 


2.00 


0.66 


186 


M52R20 


26.6 


24.0 


2975 


1910 


1.90 


0.65 


166 


M52R21 


28.8 


27.0 


3027 


1939 


2.05 


0.62 


172 


M52R22 


28.0 


29.0 


3357 


2155 


1.98 


0.61 


171 


M54R18 


24.0 


21.8 


3369 


2158 


2.07 


0.70 


224 


M54R19 


25.5 


24.8 


3556 


2279 


2.00 


0.66 


195 


M54R20 


26.3 


25.1 


3348 


2144 


2.08 


0.65 


161 


M54R21 


27.9 


25.1 


3017 


1935 


1.96 


0.63 


190 


M54R22 


30.0 


26.0 


2865 


1834 


2.15 


0.62 


198 


M60R25 


32.1 


33.0 


3699 


2367 


2.22 


0.62 


170 


Observations 




26.4±0.5 


2954±175 


2200±90 


1.27±0.44 


0.38±0.12 


198±19 



(0.64±0.15*) 



if we choose a flatter IMF for our model, it is possible to 
explain the observed slope of the mass function. 

The evolution and disruption of star clusters due to 
gas expulsion and stellar evolution may play an impor- 
tant role in the evolution of the globular cluster mass 
function (GCMF). It has been shown that these early 
evolutionary processes, i.e. early cluster dissolution due 
to gas expulsion and mass loss from stellar evolution, can 
preferentially destroy low-mass clusters and significantly 
flatten the low-mass end of the power-law initial GCMF 
( Goodwin 1997, Goodwin 1998, Kroupa & Boily 2002, 
Baumgardt et al. 2008b, Parmentier et al. 2008, Vesperini 
2009c). That is, the birth process of a globular cluster 
can be very violent, depending on the exact initial con- 
dition s of its parent gas cloud ([Baumgardt fc Kroupal 



2008 



2007; Baumgardt, Kroupa. fc Parmentier! 

Parmentier fc Gilmord 2007). Therefore, the first 



few million years can have a significant impact on 



a globular cluster's evoluti o n ov er a Hubble time. 
iMarks. Kroupa. fc Baumgardt! (|2008l ) find that, depending 
on the concentration of the initial cluster, gas expulsion 
can induce a flattening of the stellar mass function of 
the cluster. Unfortunately, this process cannot be easily 
computed in a direct way for a globular cluster of the 
pre-gas expulsive mass of Pal 14, because clusters can lose a 
large fraction of their birth stellar population as a result of 
gas expulsion and stellar evolution mass loss. We therefore 
have to use an approximation here. 

Thus, following IMarks. Kroupa. fc BaumgardtJ (|2008T ) 
and assuming that a certain flattening of the mass function 
slope has happened within the first 100 Myr of the clus- 
ter's evolution, we change the initial mass function slope 
and start the simulations at a cluster age of 100 Myr. In 
this fashion we have performed a new series of iV-body sim- 
ulations for models with various flattened initial slopes of 
the mass function instead of the canonical Kroupa IMF (see 
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Table 2. Details of all TV-body models for three sets of flattened mass functions, with initial slopes of 0.5, 0.6, and 0.7 for 
stars less massive than 0.5 M Q (see Sec. E31 for more details). Columns are the same as in Table [l] but in column 2, the 
first number in the model name is the adopted high-mass slope. Errors are the standard deviations derived from multiple 
runs with the same parameters (see text). The best fitting models that agree within the uncertainties with all observational 
parameters are highlighted with boldface. 



Model Rphm [pc] Rphl [pc] N bs M i-<S ta I M 0l a cr ;os [km/sec] N g 

(±1.6) (±2.3) (±290) (±270) (±0.13) (±0.03) (±25) 



1 


F0.5M45R18 


26.7 


24.9 


2391 


1543 


1.19 


0.59 


216 


2 


F0.5M45R20 


28.7 


27.5 


2322 


1476 


1.32 


0.56 


200 


3 


F0.5M45R22 


30.8 


29.5 


2332 


1551 


1.53 


0.53 


216 


4 


F0.5M50R18 


26.9 


27.7 


3045 


1959 


1.28 


0.62 


253 


5 


F0.5M50R20 


29.1 


27.7 


2635 


1678 


1.31 


0.60 


234 


6 


F0.5M50R22 


30.9 


30.6 


2722 


1718 


1.16 


0.57 


242 


7 


F0.5M55R18 


26.0 


26.3 


3327 


2120 


1.38 


0.66 


285 


8 


F0.5M55R20 


29.4 


27.2 


2824 


1819 


1.22 


0.62 


259 


9 


F0.5M55R22 


30.5 


29.3 


2926 


1917 


1.38 


0.60 


268 


10 


F0.6M45R18 


25.2 


25.0 


2878 


1903 


1.25 


0.61 


218 


11 


F0.6M45R20 


27.9 


25.5 


2493 


2017 


1.36 


0.59 


228 


12 


F0.6M45R22 


30.0 


28.5 


2506 


1646 


1.60 


0.56 


210 


13 


F0.6M50R18 


25.9 


23.9 


2987 


1999 


1.35 


0.64 


249 


14 


F0.6M50R20 


27.7 


28.0 


3133 


2001 


1.47 


0.61 


261 


15 


F0.6M50R22 


31.1 


31.5 


3132 


2162 


1.51 


0.60 


217 


16 


F0.6M55R18 


25.4 


23.5 


3248 


2189 


1.26 


0.68 


264 


17 


F0.6M55R20 


27.7 


26.1 


3358 


2155 


1.38 


0.64 


296 


18 


F0.6M55R22 


30.5 


29.6 


3257 


2002 


1.36 


0.63 


275 


19 


F0.7M45R18 


24.3 


21.0 


2734 


1759 


1.43 


0.68 


238 


20 


F0.7M45R20 


26.6 


25.6 


2971 


1926 


1.47 


0.61 


232 


21 


F0.7M45R22 


29.2 


29.6 


3021 


1937 


1.74 


0.59 


222 


22 


F0.7M50R18 


24.4 


22.1 


3204 


2101 


1.25 


0.68 


279 


23 


F0.7M50R20 


26.5 


25.1 


3324 


2172 


1.77 


0.65 


283 


24 


F0.7M50R22 


28.8 


28.0 


3419 


2202 


1.47 


0.69 


272 


25 


F0.7M55R18 


24.4 


21.5 


3516 


2209 


1.40 


0.72 


330 


26 


F0.7M55R20 


26.6 


25.5 


3773 


2526 


1.41 


0.68 


307 


27 


F0.7M55R22 


29.6 


27.5 


3451 


2290 


1.45 


0.66 


289 




Observations 




26.4±0.5 


2954±175 


2200±90 


1.27±0.44 


0.38±0.12 


198±19 



(0.64±0.15*) 



overview in Table [2]) . We have chosen three sets of slopes 
for the mass function: a a — {1.7,0.7}, — {1.6,0.6}, and 
a c — {1.5,0.5}, where the first number in each set is the 
slope of the mass function for stars more massive than 0.5 
Mq and the second one is for low-mass stars (for compari- 
son, the Kroupa IMF is a — {2.3, 1.3}). The second column 
in Table [2] shows the name of the simulated models. For ex- 
ample, F0.6M50R20 represents a flattened model with the 
slope of Qb, initial mass of 50000 Mg, and initial half-mass 
radius of 20 pc. 

The maximum mass in the mass spectrum was set to 
5 Mq, instead of 100 Mq. The reason can be attributed 
to the short lifetimes of the most massive stars, i.e., the 
stars with masses larger than 5 Mq will have died or turned 
into compact remnants within the first 100 Myrs. As stated 
earlier, it is reasonable to assume a retention fraction of 0% 
for compact remnants due to the low escape velocity from 
Pal 14. 

In Figure [6] we plot the mass function evolution for one 
of the best fitting models, F0.5M50R18. We see that mass 
segregation has taken place, in addition to our artificial flat- 
tening of the mass function, by about the same amount as 
in the non-fiattened case in Fig. [5] 



The results are summarized in Table [2] We find that 
a number of models in Table [5] can reproduce the present 
day properties of Pal 14. From Fig. [7| we see that all mod- 
els with a flattened IMF can reproduce the observed slope 
of the mass function much better. However, the half-light 
radius, due to the small value of the observational error, 
can constrain the calculated models. The number of bright 
stars is also well reproduced in some models. As can be seen 
in Fig. [7] models with a lower mass and slope pair a a , as 
well as models with a higher mass and slope pair a c are not 
compatible with the observed value of Nb s . 

In order to calculate the uncertainties of our results 
we have simulated each model for several times (i.e., 2 to 
6 runs for each model). This also gives us some insight to 
the uncertainties of the models in Section \4. II The resulting 
errors are also given in Table 2. These are about the same 
as the errors of the models in Table 1. 

Like in Section |4. II the line-of-sight velocity dispersion 
in most of our models is about 0.6 kms -1 and is compatible 
with the observed velocity dispersion including or excluding 
star 15. 

As can be seen in Table [2] due to the larger number 
of massive stars in the clusters with initially flattened mass- 
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Figure 4. The cumulative distribution function (cdf ) versus line- 
of-sight velocity for the sample with star 15 (upper panel) and 
without star 15 (bottom panel). The solid lines represent the 
observed data, whereas the dotted line is the cdf of one of the 
modelled clusters (M50R20). The maximum difference between 
the two cdfs is 0.13 for the upper panel and 0.14 for the bottom 
panel, which corresponds to P-values of 0.94 and 0.91, respec- 
tively. 



function, the average values of N g are larger than for clusters 
with regular initial condition. However, the number of giant 
stars in some of the best fitting models (e.g. F0.5M50R20, 
F0.6M45R18, F0.6M50R18, and F0.7M45R20) are still com- 
patible with the observed value within the uncertainties. 
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Figure 5. Mass function of main-sequence stars for one of the 
modelled clusters, 'M50R20', inside (red pluses) and outside 
(green crosses) the half-light radius, R^. The top panel shows 
the unsegregated mass function at the start of the simulation 
which was chosen to be the canonical Kroupa IMF. The observed 
mass function slope inside the half-light radius is indicated by the 
solid line. The dotted vertical lines indicate the observed range 
of stellar mass. The bottom panel shows the mass function after 
an evolution of 11 Gyr. The slope of the mass function in the 
observed range is 1.97 ± 0.15 which is larger than the observed 
value. The flattened slope of the low-mass part and the steeper 
slope of the high-mass part within the half-light radius indicate 
that mass segregation has happened in the cluster. 



4.3 Primordial mass segregation 

Another possibility to explain the fact that a standard IMF 
gives a poor fit to the observed mass function of Pal 14 
is the assumption that Pal 14 was primordially mass seg- 
regated. To test this hypothesis, we set up clusters with 
various degrees of mass segregation. The code McLuster 
allows to specify a degree of mass segregation parameter 
(hereafter, S). S = means no segregation and 5=1 
refers to full segregation, where the most massive star is 
in the lowest point of the cluster potential and the least 
massive star is in the highest point of the cluster potential 
llSubr et alJl2008h. McLuster uses the routine described in 
iBaumgardt, De Marchi, fc Kroupal l|2008T ) to segregate the 
clusters. This routine allows to maintain the desired density 
profile when increasing the degree of mass segregation while 
also asserting the cluster to be in virial equilibrium. 

For this part of the investigation we take one of our 



modelled clusters of Sec. 14.11 and add different degrees of 
primordial mass segregation. Thus, the initial 3D half-mass 
radius and initial cluster mass are fixed to M = 5OOOOM0 
and Rhm = 20 pc, respectively. The mass segregation pa- 
rameter was changed in the range S = 0.5 — 0.95. 

The results of the simulated models are shown in Ta- 
ble [3] In all cases the mass function flattens as the cluster 
loses stars. The amount of mass function flattening depends 
on the given amount of mass segregation, i.e. initially more 
segregated clusters have more flattened mass functions than 
unsegregated clusters. From Table [3] we see that low values 
of S cannot reproduce the observed mass function slope, 
while highly segregated clusters (S ^ 0.9) are able to fit the 
observed slope of the mass function. 

We also find that primordial mass segregation in clus- 
ters leads to a stronger expansion than for unsegregated clus- 
ters. Figure [8] shows the evolution of 3D half-mass radius of 
clusters for different values of the mass segregation parame- 
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Table 3. Details of all ./V-body models with primordial mass segregation and of the one model without primordial 
mass segregation but with primordial binaries. The mass segregation parameter changes in the range S = 0.5 — 0.95; 
the binary fraction for the binary model is 0.05. Columns are the same as in Table [T] just in column 1 the first 
number within the name of the model is here the adopted mass segregation parameter, S, or the fraction of stars in 
binaries, respectively. The best fitting models that agree within the uncertainties with all observational parameters 
are highlighted with boldface. 



Model 


Rphm [P c ] 


Rphi [pc] 


N bs 


M r< Rh J M el 


a 


a los [km/sec] 






±1.2 


±1.9 


±120 


±80 


±0.10 


±0.05 


±8 


S0.50M50R20 


29.4 


25.3 


3184 


2044 


1.90 


0.62 


193 


S0.60M50R20 


32.2 


27.5 


3209 


2061 


1.9 


0.58 


207 


S0.70M50R20 


32.9 


27.9 


3267 


2106 


2.0 


0.60 


191 


S0.80M50R20 


35.6 


29.6 


3219 


2066 


1.95 


0.57 


222 


S0.90M50R20 


42.4 


36.0 


3257 


2100 


1.6 


0.51 


211 


S0.91M50R20 


43.8 


38.5 


3365 


2169 


1.63 


0.46 


206 


S0.93M50R20 


47.8 


40.6 


2887 


1870 


1.27 


0.44 


217 


S0.95M50R20 


48.7 


41.5 


3159 


2045 


1.28 


0.47 


215 


S0.95M50R15 


37.0 


27.6 


3077 


1978 


1.33 


0.58 


194 


S0.90M50R15 


33.9 


26.3 


3236 


2085 


1.67 


0.55 


210 


B0.04M50R20 


27.4 


25.3 


3132 


2009 


2.0 


0.87 


185 


Observations 




26.4±0.5 


2954±175 


2200±90 


1.27±0.44 


0.38±0.12 


198±19 



(0.64±0.15*) 
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Figure 6. The same as Figure \E\ but here we started with a 
flattened mass function with a slope a = 1.5 above 0.5Mq, and 
a = 0.5 below 0.5Mq. The initial mass of this particular model is 
50000 Mq and the initial half-light radius is 18 pc. The maximum 
mass in the mass spectrum was set to 5 Mq, instead of 100 Mq. 
Further properties of the cluster after 11 Gyr of evolution are 
given in Table [2] 



ter. We see that by increasing primordial mass segregation, 
the final half-mass radius rises. For example, for a mass seg- 
regation degree varying between 0.5 to 0.95 and an initial 
half-mass radius Rh m ~ 20 pc, the final 3D half-mass ra- 
dius changes from 39 pc to 67 pc. This expansion is due to 
both dynamical and stellar evolution. The stellar evolution 
in the first 100 Myr has a more important effect and leads 
to a jump in the cluster's radius such that the cluster with 
a higher degree of segregation experiences a larger jump. 
Consequently, this increase of radius leads to a lower value 
of the mean velocity dispersion. 

In order to reproduce the final Rphi compatible with the 
observed value, we therefore have to choose a value for the 
initial 3D half-mass radius as small as Rhm = 15 pc for a 
mass segregation parameter of S ^ 0.9. In this case, the fi- 
nal line-of-sight velocity dispersion in this model is about 0.6 
kms~l, which is compatible with the observed value. Figure 
[9]shows the evolution of the mass function of this particular 
model. The mass function slope after 11 Gyr of evolution re- 
produces the observations reasonably well. This model can 
also reproduce all other observational parameters well, but 
such a high value for the mass segregation parameter im- 
plies a very tight correlation of a star's mass and its specific 
energy within the cluster at birth. There is no conclusive 
evidence for such a high degree of primordial mass segrega- 
tion in any observed star cluster. In addition, as can be seen 
in Fig. [9] the difference between the slope of the mass func- 
tion inside and outside the half-light radius is larger than 
for other models. Thus, further analysis to determine the 
observed mass function outside the half-light radius would 
help to distinguish the mass segregation scenario from other 
scenarios. 

Figure [TOl compares the projected half-light radii, R p hi, 
and the projected half-mass radii, R p hm for the three dif- 
ferent sets of modelled clusters from Sec. I4.H4.3I As can be 
seen, in almost all models R p hi is smaller than Rphm- For 
initially segregated clusters this difference is more than 10% 
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Figure 7. The half-light radius (top panel), the number of bright 
stars (middle panel), and the present day slope of the mass func- 
tion for stars with masses between 0.525 and O.795M0 (lower 
panel) for all models with flattened IMFs (Sec. I4.2j l. The hori- 
zontal axis is the model number listed in Table [2] The horizontal 
solid and dashed lines indicate the observed values and the cor- 
responding errors, respectively. The horizontal long-dashed lines 
show the total errors including the typical numerical errors given 
in Table[2] i.e. points outside these lines are not compatible with 
the observations. Further properties of the clusters after 11 Gyr 
of evolution are given in Table |21 



such that the assumption that mass traces light within a 
cluster is disputable. It must be mentioned that 10% differ- 
ence in cluster radius leads to 30% difference in the number 
of stars, which is significant and should be taken into ac- 
count in discussions. 



4.4 Primordial binaries 

Observations and theoretical work suggest that stars may 
be born with initially very large binary fractions and thus 
star clusters may contain a significant population of primor- 
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Figure 8. Evolution of the 3D half-mass radius for different 
degrees of primordial mass segregation with a canonical Kroupa 
IMF. The initial mass is 50000 Mq, and the initial half- mass ra- 
dius is 20 pc for all models, except for one model with -R^m = 
15pc. Due to stellar evolution clusters experience a rapid expan- 
sion such that clusters with higher degrees of primordial mass seg- 
regation experience a larger jump in the half-mass radius within 
the first 100 Myr of evolution. For comparison, the stellar evolu- 
tion was artificially switched off for one model. As can be seen, 
this cluster has no sudden expansion in the beginning of its evo- 
lution. 



dial binaries l|Kroupalll99"5al lbh. Low density stellar systems 
(i.e., less than 10 2 Mqpc -3 ) have a binary fraction of at 



least 50%, while in dense clusters (i.e., more than 10 



10 a 



M 



), the binary fraction is less than 10-20% l|Hut et al.l 



119921) . 

The stellar density in Pal 14 is extremely low (0.1-0.2 
Mqpc -3 ), and thus it is expected that most of the primordial 
binaries survive during the evolution of the cluster. Thus, 
the present day binary population could be close to the pri- 
mordi al binary population which it had at birth ()Kroupal 
ll995al lbt lKupper fc Kroupall2010l ) making Pal 14 a good tar- 
get to determine the primordial binary fraction of a globular 
cluster. Here we discuss one model with primordial binaries. 

Since tight binaries can significantly speed up the re- 
laxation process, one might be able to explain the flattening 
of the observed mass function, since by adding binaries, one 
increases the mean mass of particles in the cluster (as a 
binary acts like one heavy particle) , and by increasing the 
mean mass we decrease the relaxation time (|Heggie fc Hutl 
120031 ). 

Unfortunately, binaries slow down direct iV-body com- 
putations immensely because time steps have to be very 
small for their integration, such that most numerical inves- 
tigations neglect binaries completely. It is only the peculiar 
configuration of Pal 14 with its relatively low mass and low 
density that enables us to perform such a simulation with 
a few thousand binaries over a time-span of 11 Gyr. Yet, 
such a computations takes about one month GPU time on 
a GeForce GTX 280 GPU. 

We set up one of the modelled clusters of Sec. I4.1l with 
an initial mass of M=50000 Mq, initial half- mass radius of 
Rhm = 20pc, a canonical Kroupa IMF without primordial 
mass segregation but with a non-zero primordial binary frac- 



12 Hasani Zonoozi et al. 



Initial Mass Function 



6.0 
5.0 
-g4.0 

1 3.0 

2 1.0 



0.0 
-1.0 





I I 


r <^hm ' + 
r > K hm 














I 




•* ++ 

1 1 1 



-1.5 -1 -0.5 0.5 1 
log (mass) [M sun ] 

Mass Function at T=1 1 Gyr 



1.5 



5.5 - 



o 4.0 



3.5 - 



— i 1 r 

observation - 
r<R h , 




-1.2 



-0.8 -0.6 -0.4 -0.2 
log (mass) [M sun ] 



Figure 9. The same as Figure[5] but here we started with a pri- 
mordially mass segregated cluster with a canonical Kroupa IMF. 
The initial mass of this particular model (S0.95M50R15) is 50000 
Mq, the initial half-light radius is 15 pc, and the mass segregation 
parameter is set to S = 0.95. Clusters with such a strong degree 
of primordial mass segregation are able to reproduce the observed 
flat mass function inside the half-light radius (solid line). 
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Figure 10. The projected half-light radius (R p hi) plotted 
against the projected half-mass radius (R p hm) f° r the different 
sets of models listed in Tables 1-3. The black dotted line shows 
the region where Rp^i = R v hm-, whereas the cyan dash-dotted 
line demonstrates if Rp^i is 10% smaller than R p hm- Initially 
segregated clusters show more than 10% difference between R v hl 
and R p hm- 
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Figure 11. The period distribution of bound binary systems 
found in the cluster at the beginning and after 11 Gyr evolu- 
tion. Dashed and solid lines are the best Gaussian fits to the ini- 
tial and final period distributions, respectively. The coefficients 
of the Gaussian fits with a functional form of aexp[(a; — b) 2 /2c 2 ] 
are given in the figure. The errors give the 1 a errors. The pe- 
riod distribution stays the same within the errors, showing that 
dynamical encounters are unimportant in changing the binary 
population in Pal 14. Almost 85% of primordial binaries have 
survived the cluster evolution. 



tion. We add a few binaries to see if the present day binary 
distribution is a good match for the primordial one. A pri- 
mordial binary fraction, 



fbin 



Nbi; 



{N bin + N s 



(5) 



of 0.05 (i.e. about 10% of stars being in binaries because 
N s S> Nbin), where Nbin is the number of binary systems in 
the cluster and N„ is the number of single stars, is adopted, 
in order to see how many of them will be destroyed over the 
course of time and which effect on the cluster's evolution 
they have. The initial period distribution (i.e. the post-gas- 
expulsion period distribution) used in this study is the birth 
period distribution of Kroupa (1995b, 2008). Note that the 
birth period distribution is not a Gauss function. But the pe- 
riod distribution function of bound binaries (i.e., the widest 
binary systems are excluded) that results from this birth 
distribution is approximately a Gauss distribution already 
at T = Gyr because the widest initial binaries are split due 
to crowding even in the extremely low density of our initial 
Pal 14 model. 

Our simulation shows that such a small binary fraction 
does not increase the mass segregation enough, as 5% bi- 
naries do not increase the mean mass significantly. Hence, 
the obtained slope of the mass function, which is about 
a ~ 2.0 is still not compatible with the observed value, 
a = 1.27 ± 0.44 within a 1-a level of confidence. 

After 11 Gyr evolution of the cluster, the binary fraction 
reaches ft, in ~ 0.04, which means that during the evolution 
most of the binaries (~85%) survived. In other word, the 
present day binary fraction that we see in Pal 14 is almost 
equal to the primordial binary fraction at its birth. Figure 
1111 shows the period distribution of binary systems at the 
beginning and after 11 Gyr of evolution. The parameters 
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of the best-fitting Gaussian distributions to the initial and 
final period distributions are presented in Fig. If f I Since the 
mean period does not change with time (within the error), 
most of the missing binaries might have escaped the cluster 
rather than gotten disrupted. 

Moreover, adding the binaries increases the mean veloc- 
ity dispersion such that it reaches 0.87 km/sec if we consider 
all giant stars and main sequence stars with mass higher 
than 0.8 Mq. We measure 0.64 km/sec, if we add a cut- 
off in velocities leaving out the stars that are 5 km/sec away 
from the mean and hence would not be considered as cluster 
members in a cluster like Pal f 4. ft should be noted that the 
5 km/sec limit is only applied for this test to compare the ve- 
locity dispersion of our model clusters with the observations 
since observers would not count stars with such discrepant 
velocities as members. In all other cases stars were consid- 
ered as members only based on their distance to the cluster 
center. An overview of the results for this cluster can be 
found in Table [3] In order to see how well these velocity dis- 
persions are compatible with the observed sample including 
star 15, we use a KS test again. Here, the P-value is 0.68 
for both cases with and without adding the cut-off for high- 
velocity stars. This means that the data presented in Jordi 
et al. (2009) is compatible with our modelled cluster with 
primordial binaries with 68% confidence. This is compatible 
with recent results bv lKiipper fc Kroupal (|2010h . 



5 CONCLUSIONS 

We performed a series of iV-body computations in order to 
create a realistic model for the outer halo globular cluster 
Pal 14 and show that it is possible nowadays to directly 
calculate the evolution of small Milky Way globular clusters 
over a Hubble time by direct N-body simulations. For the 
sake of comparison, and as a rough estimate, the modelling 
of a typical halo globular cluster (i.e. with the initial mass of 
about 2 x 10 s Mq and half-light radius of about 5 pc) over 
12 Gyr will take about 10 times longer than Pal 14, given 
the same resources. Moreover, about half of the Galactic 
halo globular clusters are located at Rq < 6 kpc, hence 
the tidal field and cluster evaporation would be stronger, 
thereby implying a higher initial cluster mass and longer 
time of simulation. 

Pal 14's large radius together with its low-mass, thus 
low-density allow us to simulate it on a star-by-star basis. 
In this way we are not only able to learn more about the 
dynamical history of such kind of globular clusters, but are 
also able to study the process of cluster formation as this 
period leaves its imprint on its ensuing evolution. Moreover, 
we show that there is valuable synergy between simulations 
and observations, as simulations can make predictions for 
observations and, on the other hand, can test observational 
parameters on their consistency. 

Therefore, we generated 66 models divided into four 
categories: clusters with a canonical Kroupa IMF (Sec. [471]), 
clusters with a flattened IMF (Sec. 14. 2[) . clusters with 
a Kroupa IMF but with primordial mass segregation 
(Sec. 14. 3D . and finally one cluster with primordial binaries 
(Sec. 14.41) . By varying the initial half- mass radius and the 
initial mass of the model clusters, we searched through pa- 



rameter space to obtain the model that matches the avail- 
able observations of Pal 14 best. 

We used the most recent observational data to con- 
strain the modelled clusters. Thus, for all simulated clus- 
ters we compared the half-light radius, the total number of 
bright stars inside the half-light radius, and the present day 
mass function slope inside the half-light radius in the mass 
range 0.5 — O.8M0 with the observed values. Furthermore, 
we checked the line-of-sight velocity dispersion and the num- 
ber of giant stars for their consistency with the numerical 
results. 

From the computations we find that all clusters seg- 
regate significantly over the course of time, although the 
relaxation time of Pal 14 is longer than a Hubble time. We 
furthermore find that this mass segregation leads to an in- 
creasing flattening of the mass function inside half-light ra- 
dius, which is what we were looking for since the observa- 
tional value for the mass function slope of Pal 14 in the 
respective regime is only a = 1.27 ± 0.44 (|jordi et al.ll2009T ) 
compared to a Kroupa IMF value of 2.3. The initially unseg- 
regated clusters with a canonical Kroupa IMF do not pro- 
duce enough depletion in the slope of the mass function in 
order to be in good agreement with the observations, though. 
Only a few models of this set marginally agree within the 
uncertainties with this value for a. Besides that, the other 
observational constraints are matched reasonably well by 
most of the models in this set (Table [T] and Fig. [3]). 

The models with the flattened initial mass functions can 
make up for the discrepancy in a, though some of them fail 
in the other parameters (Table [2] and Fig. [7|. Primordial 
mass segregation also increases the flattening of the mass 
function but the necessary degree of primordial mass segre- 
gation turns out to be very high (Table l4~3l) . 

The initially flattened mass function can be explained 
by loss of low-mass stars from an initially mass-segregated 
cluster across the tidal radi us during the gas-expulsion 
phase (|Marks fc KroupalbOloT l. Therefore, the flattened ini- 
tial mass function of Pal 14 we uncover here strengthens 
the idea that this cluster formed in a stronger tidal field 
rather than in isolation. This may have been the case if ei- 
ther the Galactic tidal field has since then evolved, or if the 
cluster was born closer to the Galactic centre (eccentric or- 
bit), or, more speculatively, the cluster birth site was a now 
detached/disrupted dwarf galaxy. However, if Pal 14 has al- 
ways been on the same circular orbit, it is hard to imagine 
that the tidal field would have played a role during the clus- 
ter response to gas expulsion (Parmentier & Kroupa 2010). 
Primordial binaries (fun = 0.05) do not help to reproduce 
the observational value of the mass function. 

Thus, we are able to find some scenarios that can repro- 
duce all observations reasonably well. We find that the initial 
mass of Pal 14 was about 50000 Mq and that it started with 
a half-mass radius of about 20 pc. The cluster must have 
been significantly mass-segregated with stellar mass func- 
tion depleted in low-mass st ars. These are the signatures of 
a pos t-gas expulsion cluster (jMarks. Kroupa. fc Baumgardtl 
120081 ). The inferred post-gas expulsion "initial" half-mass ra- 
dius is larger than the half-mass radii of most globular clus- 
ters which is of the order of 3 pc. Moreover, we show that the 
different scenarios show different mass functions outside the 
half-light radius, such that finding the observational value 
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for this region would help to discriminate between the dif- 
ferent scenarios. 

The effect of unresolved binaries is an important issue 
that has to be considered in the observational values of the 
velocity dispersion and the mass function. If Pal 14 has a 
significant binary population, which we expect due to the 
low density of this cluster, the true (binary-corrected) mass 
function should be steeper, because of the existence of unre- 
solved components (jKroupa. Gilmore. fe Toutlll99ll ). There- 
fore, the slope of the mass function is the most uncertain pa- 
rameter of the observations of Pal 14. Even without taking 
binaries into account Jordi et al. (2009) estimate the un- 
certainty to be ~ ±0.4. Furthermore, as Kiipper & Kroupa 
(2010) have shown, if the binary fraction was normal (i.e. 
more than 50%), the binary-corrected velocity dispersion of 
the cluster would also be smaller than the value reported by 
Jordi et al. (2009). 

From the set of computations we furthermore find that 
primordial mass segregation leads to a larger impact of stel- 
lar evolution on cluster expansion. This is, segregated clus- 
ters are more affected by mass loss through stellar evolution 
(Fig. [8]). We furthermore quantify the effect of mass segre- 
gation on the ratio of half-mass radius to half-light radius 
and find the latter to be up to 10% smaller (Tig. llO[) . leading 
to significant observational uncertainties for the number of 
stars within this radius in such cases. 

Moreover, we measured the line-of-sight velocity disper- 
sions of our models taking into account all stars, and com- 
pared it with the observed values. We additionally measured 
the velocity dispersion considering only giant stars, and find 
that both values do not show significant difference between 
each other, showing that giant stars are good tracers of a 
cluster's global velocity dispersion. Most clusters show a ve- 
locity dispersi on of ~ 0.6 k m/se c, suggesting that starl5 of 
the sample of iJordi et al. I (12009T ) might be a normal cluster 
member. 

Finally we find from the model with a primordial binary 
fraction of 5% that, due to the low stellar density of Pal 14 
the initial binary distribution stays almost unchanged by the 
dynamical evolution of the cluster, and only about 15% of 
binaries escape or get disrupted in the course of time (Fig. 
Hip . Low-density clusters like Pal 14 are therefore good ob- 
jects to study the primordial binary distribution of globular 
clusters. 



ACKNOWLEDGEMENTS 

A.H.Z and H.H thank the Stellar Population and Dynamics 
Group of the Argelander Institute for Astronomy and the 
Institute for Advanced Studies in Basic Sciences (IASBS) 
for providing financial support for this research. 



REFERENCES 

Aarseth S. J., 2003, Gravitational N-Body Simulations, 
Cambridge, UK, Cambridge University Press 

Anders, P.; Lamers, H. J. G. L. M.; & Baumgardt, H. 2009 
A&A, 502, 817. 

Baumgardt H., Makino J., 2003, MNRAS, 340, 227. 



Baumgardt H., Grebel E. K., & Kroupa P., 2005, MNRAS, 
359, LI 

Baumgardt H., Kroupa P., 2007, MNRAS, 380, 1589 
Baumgardt H., De Marchi G., Kroupa P., 2008a, ApJ, 685, 
247 

Baumgardt H., Kroupa P., Parmentier G., 2008b, MNRAS, 
384, 1231 

Bonnell I. A. & Davies M. B., 1998, MNRAS, 295, 691 
De Marchi G., Paresce F., Pulone L., 2007, ApJ, 656, L65 
Dekker H., D'Odorico S., Kaufer A., Delabre B., & Kot- 

zlowski H., 2000, Proc. SPIE, 4008, 534 
Djorgovski S. & King I. 1986, ApJ, 305, L61L65. 
Dotter A., Chaboyer B., Jevremovic D., Baron E., Ferguson 

J.W., Sarajedini A., & Anderson J., 2007, AJ, 134, 376 
Dotter A., Sarajedini A., & Yang S.-C., 2008, AJ, 136, 1407 
Elmegreen B.G., 2004, MNRAS, 354, 367 
Giersz M. & Heggie B.C., 2009, MNRAS, 395, 1173. 
Goodman J. & Hut P., 1989, Nature, 339, 40. 
Goodwin S. P., 1997, MNRAS, 286, 669 
Goodwin S. P., 1998, MNRAS, 294, 47 
Haghi H., Baumgardt H., Kroupa P., Grebel E. K., Hilker 

M., Jordi K., 2009, MNRAS, 395, 1549 
Harfst S., Portegies Zwart S., Stolte A., 2009 

(|arXiv:0911.3058|l 
Harris W. E., 1996, AJ, 112, 1487 

Heggie D. C. & Hut P., 2003, The Gravitational Million- 
Body Problem, Cambridge: Cambridge University Press. 
Hilker M., 2006, A&A, 448, 171. 

Hurley J. R., Pols O. R., Tout C. A., 2000, MNRAS 315, 
543 

Hut P., et al., 1992, PASP, 104, 981 

Jordi K., Grebel, E. K., Hilker, M., Baumgardt, H., Frank, 
M., Kroupa, P., Haghi, H., Cote, P., Djorgovski, S. G., 
2009, 137, Issue 6, pp. 4586-4596. 

King I. R., 1966, AJ, 71, 64 

Kruijssen, J. M. D.; Mieske, S., 2009 A&A, 500, 785 
Krumholz M. R., Klein R. I., McKee C. F., Offner S. S. R., 

Cunningham A. J., 2009, Science, 323, 754 
Kroupa P., Gilmore G., Tout C.A., 1991, MNRAS, 251, 293 
Kroupa P., 1995a, MNRAS, 277, 1491 
Kroupa P., 1995b, MNRAS, 277, 1507 
Kroupa P., 2008, 2008, The Cambridge N-Body Lectures, 

Lecture Notes in Physics, 760, 181 
Kroupa P., 2001, MNRAS, 322, 231 
Kroupa P., 2002, Science, 295, 82 
Kroupa P., Boily C. M., 2002, MNRAS, 336, 1188 
Kiipper A.H.W., Kroupa P., 2010, ApJ, 716, 776 
Marks M., Kroupa P., Baumgardt H., 2008, MNRAS, 386, 

2047 

Marks M., Kroupa P., 2010, MNRAS, 406, 2000 
Mouri H., Taniguchi Y., 2002, ApJ, 580, 844 
Parmentier G., Gilmore G., 2007, MNRAS, 377, 352 
Parmentier, G., Goodwin, S. P., Kroupa, P., Baumgardt, 

H., 2008, ApJ, 678, 347 
Parmentier G., Kroupa P., 2010, accepted in MNRAS 
Salpeter E. E., 1955, ApJ, 121, 161 
Sarajedini A., 1997, AJ, 113, 682 

Subr L., Kroupa P., Baumgardt H., 2008, MNRAS, 385, 
1673. 

Tan J. C, Krumholz M. R., McKee C. F., 2006, ApJ, 641, 
L121 



Numerical Modelling of Palomar 14 

Trenti M., Vesperini E. & Pasquato M., 2010, ApJ, 708, 
1598 

Vesperini E. & Heggie D. C, 1997, MNRAS, 289, 898 
Vesperini E., McMillan S., Portegies Zwart S., 2009a, 
Ap&SS, 177 

Vesperini E., McMillan S., Portegies Zwart S., 2009b, ApJ, 
698, 615 

Vesperini E., 2010, Philosophical Transactions of the Royal 
Society A: Mathematical, Physical and Engineering Sci- 
ences, 368, 829 l|arXiv:0911.0793|> 

Vogt S. S., et al., 1994, Proc. SPIE, 2198, 362 



This paper has been typeset from a TffjX/ M?gX file prepared 
by the author. 



